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We study the spatio-temporal dynamics of the generalized repressilator, a system of 
coupled repressing genes arranged in a directed ring topology, and give analytical 
conditions for the emergence of a cascade of unstable periodic orbits (UPOs) that 
lead to reachable long-lived oscillating transients. Such transients dominate the 
finite time horizon dynamics that is relevant in confined, noisy environments such as 
bacterial cells (see our previous work [Strelkowa and Barahona, 201 0J) and are there- 
fore of interest for bioengineering and synthetic biology. We show that the family of 
unstable orbits possesses spatial symmetries and can also be understood in terms of 
traveling wave solutions of kink-like topological defects. The long-lived oscillatory 
transients correspond to the propagation of quasistable two-kink configurations that 
unravel over a long time. We also assess the similarities between the generalized 
repressilator model and other unidirectionally coupled electronic systems, such as 
magnetic flux gates, which have been implemented experimentally. 
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Bacterial cells are noisy environments that evolve over a limited time horizon. 
Under these conditions, the relevant dynamics is not dictated exclusively by 
the long-term attractors of the system but can be dominated by reachable long- 
lived transients with qualitatively different observable dynamics. Such transients 
are an important feature of naturally occurring biological networks, in partic- 
ular in developmental biology, and can also be exploited for forward design of 
genetic elements in Synthetic Biology. Here, we investigate the dynamics of 
unidirectional rings of genetic repressors and show that their observable dy- 
namics is dominated by oscillatory transients around unstable periodic orbits 
with strong spatio-temporal symmetries. The oscillatory transients around such 
orbits can be seen to correspond to the propagation of groups of kink-like ex- 
citations against the background of a dimerized solution, a picture with strong 
parallelisms with classic discrete lattice models^HS, some of which have been im- 
plemented experimentally with electronic elements. 



I. INTRODUCTION 

The dynamics of symmetric interconnected systems are an intricate reflection of the in- 
herent invariance associated with the symmetry present in the systenP. Engineered systems 
often resort to the use of symmetry because of the simplification in design and construction 
they afford, but also as a means to the production of controlled outcomes. Loops or ring 
structures of nonlinear elements®^ constitute one of the simplest topologies of practical in- 
terest in biology and engineering®^. For instance, in a recent series of papers E3 9 l 12 l 13 l Bulsara 
et al reported the implementation of rings of unidirectionally coupled magnetic flux-gates 
through electronic circuits and their successful experimental use as highly sensitive magnetic 
field sensors^* 9 * 12 * 13 !. The mathematical analysis of the steady-state dynamics of a variety of 
such ring systems has also revealed their rich dynamical behavior, including multi-stability, 
unstable periodic solutions, quasi-periodicity and even chaos and hyperchao ij 7 1 9 1 13 ! 

These classic design concepts have found new applicability with the recent emergence of 
Synthetic Biology, a field which has put to use ideas from electrical engineering, computer 
science and physics to engineer circuits made up of biological parts, in many cases through 

2 



UPOs in the generalized repressilator 



the design of gene regulatory network d 14 * 15 l Such synthetic circuits of interconnected genes 
have been implemented in living microbes to perform low-level functions, such as logical 
gates^, switches^ and genetic oscillators^, exploiting the steady-state behavior of the gene 
network. 

Here, we study the spatio-temporal dynamics of the generalized repressilator, a gene 
regulatory network of repressing genes cyclically coupled in a directed ring topology (see 
Fig. [l^,), and we concentrate on the investigation of the role of unstable periodic orbits 
(UPOs) on the observable transient dynamics. The study of the generalized repressilator 
model has a long history^^. As shown analytically in [19l there is a fundamental difference 
between rings with odd and even number of monotonic repressors: at steady state, odd- 
numbered rings possess a limit cycle attractor while even-numbered rings are bistable and 
act as switches. Hence the generalized repressilator model bridges between switches and 
oscillators in a natural way. In fact, the first ever genetic oscillator realized in a living 
bacterium by Elowitz and Leibler^ is a special case of the generalized repressilator for three 
genes, while the first genetic switch implemented in E. coli by Gardner et alJ^ corresponds 
to the generalized repressilator with two genes. 

Previous mathematical work on these systems was mainly devoted to stability analysis 
and the characterization of the attractors that determine the long term behavior of the 
rings. On the other hand, the study of the transient behavior, in particular in relation to the 
presence of UPOs, was broadly ignored. However, the standard steady-state analysis alone 
is not sufficient in systems where transients are qualitatively different from steady states 
and extend over long times compared to the observational time scales of the system, e.g., 
the average life time of a bacterial cell for a system like the repressilator. Even numbered 
repressilator rings are examples of such systems: although their steady state is typically 
bistable, with two attracting fixed points^, the fixed points are approached through long- 
lived oscillating transients that dominate the observable dynamics^. In a recent paper, we 
showed that this transient behavior could be exploited to produce switchable oscillations by 
using feedback to operate the system around these unstable oscillations^. 

In this paper, we show that the observed long-lived transients have their origin in the ex- 
istence of UPOs that attract initial conditions to their vicinity, thus funneling the dynamics 
into a lower-dimensional oscillatory subspace from which it is difficult to escape. The UPOs 
appear through a cascade of Hopf bifurcations, for which we provide analytical conditions, 
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FIG. 1. The generalized repressilator model (a) The network consists of genes interconnected 
in a directed ring topology, in which the protein product of gene j represses the expression of gene 
j + 1 with periodic boundary conditions, (b) As the parameter c (see Eq. [3]) is increased, a cascade 
of Hopf bifurcations emerges — here shown through p\ vs. P2 phase plots. In odd rings (exemplified 
by N = 9), an attracting limit cycle (meshed style) emerges from the first Hopf bifurcation (H). 
In even rings (exemplified by N = 10), a pitchfork bifurcation (P) leads to two symmetrical stable 
fixed points (solid black lines), which we refer to as 'up/down' solutions. In both even and odd 
rings, further Hopf bifurcation lead to emergence of UPOs. The difference in the bifurcation 
cascades between even and odd gene rings can be understood from the pattern of the eigenvalues 
of the linearization around the symmetric fixed point p m . Some limit cycles in even rings undergo 
Neimark-Sacker bifurcations (NS) leading to an unstable torus. 



and possess spatio-temporal symmetries consistent with tilings of the ring with smaller basic 
solutions^. In addition, the dynamics around the UPOs can also be understood in terms 
of traveling kink-like excitations propagating against the backdrop of a dimerized discrete 
lattice — a picture analogue to other propagating waves in classical physics^^. The reach- 
ability and persistence of the dynamics around the UPOs, which make them observable, is 
related to the interaction and propagation of moving kinks in the lattice. Finally, we show 
that these features also apply to the analysis of rings of magnetic flux-gates^ implemented ex- 
perimentally with electronics — a system in which extremely long-lived oscillating transients 
(especially in the presence of noise) have been observed both experimentally and numerically 
in even-numbered rings that should be bistable according to steady state analysis. 
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II. BIFURCATION ANALYSIS OF THE GENERALIZED 
REPRESSILATOR AND UPOS 

A standard model of the generalized repressilator with N genes is given by the following 
ODEs: 

™>j = c i/(Pi-i) ~ c i™j 

pj = c 3 mj - apj j = 1, . . . , N. (1) 

Here, pj and rrij describe protein and mRNA concentrations for each gene (with periodic 
boundary conditions such that p = p^). The kinetic constants c\ and c 3 represent the 
strength of repression and transcription, respectively, while c<i and C4 correspond to degra- 
dation of each species. The repressing behavior is encapsulated by f(p), a positive function 
/ : 1R + — > 1R + that is bounded between and 1 and is monotonically decreasing over its 
domain. For such a function, the system (JIJ has a spatially uniform fixed pointful in which 
all proteins have the same concentration: pj = Pj~i = p m Vj. 

Throughout this paper, our numerics are performed using a standard Hill-type function, 
which is well established in biochemical modeling: 

/GO = i^e. (2) 

where h is the Hill coefficient. However, the outlined stability analysis is similar for any 
positive function that is bounded and monotonically decreasing. For the Hill function ([2]), 
the uniform fixed point solution p m follows from solving the polynomial equation: 

Pm+P^^ — ^C. (3) 
C2C4 

For the case h = 2, which is relevant in some Synthetic Biology implementations^, the 
cubic equation that ensues can be solved explicitly using the Cardano method. Note that 
this solution branch is independent of the number of genes in the ring and exists for all 
positive c. In fact, p m is the only fixed point when the repression (or coupling) c is small 
enough. 

As the strength of the repressive coupling c increases, there is a qualitative difference 
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between the bifurcations undergone by even- and odd- numbered rings (see Fig. [I]). In odd 
rings, p m becomes unstable via a Hopf bifurcation resulting in a stable limit cycle. Even 
rings, on the other hand, first undergo a pitchfork bifurcation where two branches of stable 
dimerized solutions emerge. No further stable solutions appear as c is increased either for 
even or odd rings. Hence from the viewpoint of strict stability analysis, odd rings behave as 
oscillators, with a globally attracting limit cycle, and even rings behave as switches, due to 
their bistability. 



There are however other bifurcations that lead to the emergence of unstable solutions. 
In both even and odd rings, further increasing the coupling c leads to a cascade of Hopf 
bifurcations, where each bifurcation gives rise to a UPO that emerges from the already 
unstable branch p m . In addition, numerical continuation^ a l so reveals the existence of 
Neimark-S acker bifurcations in particular UPOs of even rings, in which an unstable torus 
emerges from the UPO. Figure [T] shows characteristic examples of these behaviors. 

The conditions for the appearance of UPOs can be derived analytically from the instabil- 
ities of the spatially uniform fixed point p m . Eliminating the mRNA variables (rrij), Eq. ([I]) 
can be rewritten as: 

Pj = ciczf{pj-i) - (C2 + c 4 )pj - c 2 c 4 pj. (4) 
This corresponds to the system 

Pj = Qj 

% = CiC 3 /(p J _i) - (c 2 + c 4 )qj - c 2 c 4 pj, (5) 
with linear and nonlinear parts as follows: 

z = (I N ® A)z + H(z). (6) 

Here, z = (. . . pj qj . . .), In is the iV-dimensional identity matrix, A = I 

\ -c 2 c 4 -(c 2 + c 4 ) 

and H(z) is a vector that contains the nonlinear coupling between the elements. 
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Linearizing around the fixed point p m we obtain the Jacobian 



J\ Pm = {In®A) + {G®H*) 



(7) 



where H* = c\c^ 







dp 







and the directed ring topology of the network is encoded in 



Pm 



the N x N adjacency matrix 



G 



... 1 
10 ... 



... 1 
\0 ... 1 ) 

Fourier-diagonalizing^SI the circulant matrix G, we obtain the following complex quadratic 
equation for the eigenvalues Xg of the Jacobian ([7]): 



A| + (c 2 + Ca)\i + C 2 C 4 1 



dp 



e n 







(8) 



with £= {0,..,iV-l}. 

The necessary condition for a Hopf bifurcation, i.e., the existence of purely imaginary 
roots, leads to inequalities on the parameters^ that determine the number of bifurcations. 
It then follows that the maximum number of Hopf bifurcation is (N/2 — 1) for even rings 
and (N — l)/2 for odd rings (see Appendix). 

If the repressor is described by the Hill function ([2]), we have 



d[_ 
dp 



hp. 



h+l 



Pm 



When h = 2 and with parameters relevant for the reported Synthetic Biology experiments^, 

for odd rings and [| f y — ljj for 



1 (JV+l) 

2 2 



the number of Hopf bifurcations is limited to 
even rings (see Appendix and Table [T]) . 

Apart from the first bifurcation in odd rings, which gives rise to an attracting limit cycle, 
all other Hopf bifurcations give rise to unstable periodic orbits (UPOs), as shown in Figure [l] 
and summarized in Table |TJ Although UPOs do not lead to long term attractors, they can still 
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be important for practical applications if they are likely to be reached from random initial 
conditions (or as a result of random fluctuations) and they dictate long-lived dynamics that 
will therefore feature in finite-time observations. As shown below, the UPOs that emerge 
from the first Hopf bifurcation in even rings have these characteristics for parameter ranges 
consistent with biological experiments. In particular, with mRNA and protein life-times on 
the order of minutes and Hill coeficients h G {2, . . . ,6}, the transient oscillatory behavior 
would extend over several hours in rings of moderate size N ~ 10 when the cellular doubling 
time is on the order of one hour. Therefore the transient behavior dominates the observable 
dynamics. 

We now explore the different dynamical facets of the emerging unstable periodic orbits. 
In particular, we describe their spatio-temporal symmetry, and their relationship to tilings 
of the ring, and present their linear stability. We also show that these solutions can be 
understood in terms of traveling waves that correspond to propagating interacting kinks, 
analogous to solutions from classical physics such as domain walls. We also show that, 
although unstable, these periodic solutions are reachable from random initial conditions and 
long-lived and, as such, they are relevant in finite-time dynamics in real-world applications. 

III. SPATIO-TEMPORAL SYMMETRIES OF UPOS OF THE 
GENERALIZED REPRESSILATOR MODEL 

The discrete rotational symmetry of the equations (Jl| induces non-trivial symmetries 
in the UPOs which can be directly related to the linear stability of the solutions. Table [I] 
presents a summary of the bifurcation analysis and emerging solutions for N e [3, 18] carried 
out with the continuation software MATCONT^. The analysis reveals that some of the 
higher order bifurcations for larger rings have the same periods as lower order bifurcations 
in smaller rings. This is the result of a rotational symmetry of such solutions that correspond 
to a tiling of the ring, i.e., gene pj has exactly the same time behavior as gene Pj+ S . This can 
be seen, for instance, in the unstable periodic solution that emerges from the second Hopf 
bifurcation in the N = 12 gene ring, where pj,pj +3 ,pj +6 and pj + g have the same time traces 
(see Fig. [2] a). Each triplet of consecutive elements follows a dynamics equivalent to that of 
a 3-gene ring, i.e., the solution corresponds to the tiling of the 12-gene ring by four 3-gene 
repressilators. Therefore, its period is the same as that of the 3-gene ring for the same set 
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of constants. 
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FIG. 2. Spatio-temporal symmetries in the generalized repressilator model, (a) Phase 
plots for the symmetric UPO originating from the second Hopf bifurcation in the 12-gene ring. We 
plot the time evolution of the first gene against other genes in the ring. Genes that are dynamically 
equivalently are plotted column-wise. Each of the four rows shows a group of genes, which behave 
dynamically in the same way as the 3-gene ring. The solution has four kinks corresponding to the 
domain walls between the 3-gene tiles, as shown in the bottom diagram, (b) Phase plots for the 
non-symmetric UPO originating from the first Hopf bifurcation in the 16-gene ring. In this case 
there is no obvious symmetry between the different genes. As shown in the bottom diagram, this 
a basic solution that contains two kinks. 



Table [I] summarizes the Hopf bifurcations and some of the properties of the ensuing 
periodic solutions (limit cycles and UPOs) obtained through numerical continuation and 
from the analysis of numerical simulations. We have considered rings of length between 
N = 3 and 18 of the model with Hill-type repression with h = 2 for a range of the bifurcation 
parameter c\ G (0,3]. The observed bifurcations follows the number and pattern of Hopf 
bifurcations found analytically (see Appendix). The solutions are sorted by increasing period 
and therefore symmetric tiling solutions appear grouped together (since they have the same 
period). Not all the solutions correspond to tilings of the ring. Some of the UPOs are 
'basic solutions' with no strict rotational symmetry, i.e., no two genes behave in the same 
way (Figure^). The fact that some solutions of the system can be understood in terms of 
tilings of basic solutions indicates the importance of domain walls and kink-like excitations in 
the understanding of the observed solutions. We explore this alternative dynamical picture 
below. 
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N 



11 

18 

■ 

17 

10 

13 

16 

3 

6 

9 

12* 
15 
18 
17 
14 
11 

£ 

13 
18 



17 
12 

I 

16* 
9 

18 
11 
13 
15 
17 



Unstable 
Floquet 



I 



I 



4/22 
7/36 
2/14 
5/28 

6/34 
3/20 
4/26 
5/32 
0/6 
1/12 
2/18 
3/24 
4/30 
5/36 
4/34 
3/28 
2/22 
1/16 
3/32 
2/26 
3/36 
0/10 



MO/10 
1/20 
2/30 



2/34 
1/24 
0/14 
1/28 

1/32 
0/18 
1/36 
0/22 
0/26 
0/30 
0/34 



Period 



86 

89 

94 

94 

100 

104 

110 

113 

132 

132 

132 

132 

132 

132 

154 

160 

169 

186 

186 

204 

213 

239 

239 

239 



274 
291 
342 
342 

393 
444 
444 
545 
645 
746 
847 



No. of 
kinks 



6 

7 

4 
5 
6 
1 

2 
3 
4 
5 
6 
5 
4 
3 




2 

2 
1 
2 
1 
1 
1 
1 



Basic solution 
(kinks) 



11(5) 
18(8) 



7(3) 



2 8(2) 
4 8(2) 



17(7) 
10(4) 
13(5) 
16(6) 
3(1) 
3(1) 
3(1) 
3(1) 
3(1) 
3(1) 
17(5) 
14(4) 
11(3) 
8(2) 




7(1) 

16(2) 

9(1) 

9(1) 

11(1) 

13(1) 

15(1) 

17(1) 



No. of 
tilings 






TABLE I. Bifurcations, symmetries and kinks in the generalized repressilator model. 

The analytical conditions (Eq. [8]) for the emergence of UPOs are confirmed for rings of length 
N = {3, . . . , 18} through simulations and numerical continuation^!. The numerics have been 
carried out for the Bill function with h = 2 for the following set of parameters: C2 = 0.12 min -1 , C3 = 
0.16 min -1 , C4 = 0.06 min -1 , and the bifurcation parameter c\ in the range (0, 3]. We have followed 
the UPOs numerically calculating their Floquet exponents, period, and spatio-temporal symmetry. 
The number of kinks has been numerically extracted from time traces according to Eq. Q. The 
UPOs are sorted by increasing period and solutions whose dynamics are equivalent by tiling are 
highlighted with the same color. For instance, the basic 3-gene repressilator (bold red text) has an 
associated family of symmetric tiling solutions for rings where the gene number N is a multiple of 
three. Phase plots of the two solutions marked with asterisks are presented in Figure [2j 
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IV. UPOS AS TRAVELING KINK SOLUTIONS 

An alternative viewpoint to understand the dynamics of gene rings is to view the periodic 
solutions in terms of propagating kinks against the backdrop of a dimerized ('up-down') 
solution. Analogues to this picture are abundant in both the nonlinear waves literature and 
in classical examples in condensed matter physics, such as propagation of domain walls in 
antiferromagnetic chains or soliton-like excitations in dimerized polymer chains^^. 

The link between the dynamical bifurcations presented above and the kink-domain pic- 
ture follows straightforwardly from the examination of the fixed points of the system. As 
explained above, the uniform solution, where all proteins have the same value pj = p m , Vj, 
is always a fixed point of the system and is the global attractor of the equations for small 
enough values of the coupling parameter c. 

When the coupling between the genes is large enough, the uniform solution becomes 
unstable either through a Hopf bifurcation (for odd rings) or through a pitchfork bifurcation 
(for even rings). The two stable fixed points that emerge from the pitchfork bifurcation 
in even rings correspond to dimerized solutions ('up-down') where genes are over-expressed 
(p u >Pm) or under-expressed (pd<p m ) in an alternating sequence around the ring. This 
arrangement leads to two degenerate solutions, equivalent by a one-gene rotation, which 
correspond to the sequences 'up-down' and 'down- up', respectively. Clearly, such a dimerized 
solution is compatible with an even ring but not with an odd ring. 

In the case of an odd ring, trying to fit an 'up-down' solution leads to frustration due 
to the boundary conditions, originating a kink mismatch when wrapping up the sequence. 
Indeed, the first Hopf bifurcation in odd rings can be viewed as the emergence of a prop- 
agating solution of this type. Hence the attracting limit cycles in odd rings correspond to 
the fundamental solution with one traveling kink propagating against a background of the 
dimerized solution. The localized structure of the kink (domain wall) is shown in Figure |3^i 
for a ring of size N = 51. A quantitative description of the dynamics of the kinks can be 
obtained from numerical simulations by defining the threshold function: 



which allows us to define the domain regions as sequences of alternating genes. The location 




1 if pj < p., 



if Pj >Pm = 'up' 



fit 



= 'down', 



(9) 
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of the kinks is then defined by the condition v(pj) = v(p 



3+1 J 



P, 




FIG. 3. Kinks in the repressilator model. The kink-domain picture provides an alternative 
characterization of periodic orbits. Numerical time snapshots in gene space (black dots) are mapped 
onto domains and kinks using Eq. [9] (gray arrows), (a) Frustrated solution in an odd- numbered 
(N=51) gene ring, which emerges from the first Hopf bifurcation and is a stable attractive limit 
cycle. Its stability follows intuitively from the fact that the kink can propagate along the ring 
forever, (b) The solution with two kinks in the even-numbered ring (N=50) is, on the other hand, 
unstable. This solution originates from the first Hopf bifurcation in an even ring. The dynamics 
around this UPO corresponds to the motion of the two kinks along the ring with a slightly different 
speed. As one kink gains on the other, they eventually annihilate and the ring becomes quiescent. 
Importantly for applications, this happens only after a very large number of periods (~ 100). 

The pattern now becomes clear. The cascade of Hopf bifurcations corresponds to the 
emergence of solutions with an increasing number of kinks in the ring. The first Hopf 
bifurcation in even rings corresponds to the emergence of a solution with two traveling 
kinks which can be seen as localized excitations in the ring, as seen in Figure [3]d for a 
ring of size N = 50. The kinks (or walls) originate from the mis- matching of two 'up-down' 
domains to give an arrangement that can be fitted into an even ring, as shown at the bottom 
of Figure [2}d. As indicated by Floquet analysis, this solution is an unstable periodic orbit 
in which the two kinks get closer as time evolves until they annihilate each other. As the 
two kinks conflate, the domain separation unravels and we recover the 'up-down' dimerized 
solution. Similarly, the second Hopf in odd rings gives rise to UPOs that correspond to 
solutions containing three traveling kinks, while the second Hopf in even rings leads to the 
emergence of UPOs with four propagating kinks in the ring. 

It then follows that even numbered rings only admit (unstable) configurations with an 
even number of kinks that eventually annihilate to recover the stable dimerized solution. 
In odd rings, only solutions with an odd number of kinks are possible and, similarly, all 
the kinks but one eventually annihilate, leaving a single kink propagating along the ring 
indefinitely, i.e., we recover the stable limit cycle basic solution. Note also that when the 
solutions correspond to tilings of basic solutions, the total number of kinks is the result of 
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multiplying the number of kinks of the basic solution by the number of tilings. Moreover, the 
number of Hopf bifurcations, which follows from symmetry constraints, is therefore linked 
to the maximum number of kinks that can be accommodated in the rings. A summary of 
the periodic solutions in Table |T] is presented in terms of their kink density in Figure |4j 

1/N 2/N 3/N 4/N 5/N 6/N 7/N 8/N 

5 



10 




9 7 5 4 3 7 

kink density 

FIG. 4. Tilings and kink density of UPOs of the generalized repressilator model. 

Summary of the solutions presented in Table [I] with the same parameters. For this set of parameters 
the kink density is always smaller than 1 /2. Solutions with the same number of kinks are connected 
by the grey (odd) and pink (even) curves. The solutions of odd-numbered rings appear at the 
crossings of the grey horizontal lines with the grey curves. Similarly, the intersections of the pink 
horizontal lines and the pink curves mark the UPOs for even rings. The colored vertical lines 
connect the families of symmetry tilings (colors matched to those in Table [I]) . 



V. OBSERVABLE LONG-LIVED OSCILLATING TRANSIENTS AROUND 
UPOS 

The kink-domain picture not only provides insight into the existence and structure of 
periodic solutions, but also gives an intuitive understanding of their stability. As expected, 
the larger the number of kinks, the more unstable the solutions become. This can be seen 
in the fact that the number of unstable Floquet directions grows as the number of kinks 
increases: the number of unstable Floquet directions is equal to the number of kinks minus 
one (Table [T]). In addition, the increased instability of the solutions with a large number of 
kinks also translates into a reduced length of the transients observed numerically, i.e., the 
time it takes for the kinks to conflate decreases as the density of kinks increases. This means 
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that for practical purposes, the dominant long-lived oscillatory transients will be associated 
with the UPO emerging from the first Hopf bifurcation in even rings. We show below that 
the length of the oscillatory transients is long compared to the observation time scales of 
the system, hence these UPOs dictate the observable dynamics of the system. 

Slow convergence to the dimerized solution through kink-kink annihilation 

Even though even rings act as bistable switches at steady state, the UPOs emanating 
from the first Hopf in even rings have a strong influence on their observable dynamics. Often, 
when starting from any random initial condition, even rings of moderate length get trapped 
into long-lived oscillating transients before they eventually approach the stable 'up-down' 
solution. These observed transients in the numerics can be traced to dynamics around the 
first unstable UPO, a solution with two kinks in the ring. A typical transient will quickly 
settle down into a configuration with two kinks separated by a distance that depends on 
the initial condition. The two kinks then propagate along the ring for a large number of 
periods, slowly approaching each other, until they annihilate. The dynamics of this long 
transient is shown in detail in Figure [5] from different perspectives. The very slow slip of 
the kinks is observable not only through the persistence of the inter-kink distance during 
the transient (Fig. J5Jd,cI) but also in the phase plot and Poincare sections of the dynamics 
(Fig. (5}^,e). In particular, the structure of the phase plot (Fig. [5^) indicates the likelihood 
that the dynamics of the transient is restricted to a reduced manifold, possibly the unstable 
torus that emerges from a Neimark-Sacker bifurcation (see Fig. Wp)- 




Observability of the transient dynamics around UPOs: reachability and 
persistence 

For the transient dynamics described above to be relevant for applications and design 
it needs to be observable, i.e., it should be reachable from large regions of phase space 
and it should persist over long time scales as compared to the intrinsic time horizon of the 
system. Our numerical simulations indicate that the UPOs that emerge from the first Hopf 
bifurcation in even rings (corresponding to the two-kink arrangement) possess these two 
attributes. Figure [6] shows that the probability of observing long oscillating transients is 
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FIG. 5. Slow kink-kink annihilation leads to long oscillatory transients. Dynamics around 
the UPO arising from the first Hopf bifurcation in the 16-gene ring to illustrate the long oscillating 
transient that eventually reaches the stable 'up-down' solution, (a) Starting from a random initial 
condition, the ring quickly settles into a configuration with two kinks, (b) The distance between the 
kinks, kdi goes from = 5 to 4 and quickly through 3, 2 until the kink-kink annihilation producing 
a long oscillatory transient of ~ 80 periods, (c) Phase plot {p\ versus p§) from the initial condition 
(light green) to the final point at the stable up-down solution (red). The black dots correspond to 
Poincare crossings at equidistant time intervals (every period) and illustrate the very slow initial 
drift which only accelerates rapidly at the end of the transient, when the unravelling of the solution 
takes place. Note also the uniform spreading of the trajectories on the surface of the phase plot, 
an indication of the reduced dimensionality of the transient dynamics, (d). Space-time diagram 
of the transient shows the slow approach of the kinks and their annihilation. The kink traces are 
indicated by gray lines, (e) Poincare sections in gene space. The value of each protein at the 
Poincare section is represented as a phasor expi2n(j)j where <j>j = Pj/p-max- The kink annihilation 
event in this representation corresponds to the unravelling of a "twist", which is a global event 
involving several ring elements and therefore occurring on very slow time scales. 
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higher than 50% for moderately long rings. In fact, all the oscillatory transients reached 
from random initial conditions correspond to the two-kink traveling solution, so this is the 
only UPO that is observable. Note also that, as expected, the longer the ring, the more 
likely it is to observe these transients since it is easier for the initial condition to settle onto a 
configuration with two well-separated kinks in the ring. This is also seen in Figure [7], where 
we show the increasing length of the oscillatory transient from a similar initial condition as 
we increase the length of the ring. These observations explain why although UPOs with 
more than two kinks are possible in even rings, they are not observed numerically when 
starting from random initial conditions. Multi-kink solutions are both combinatorially less 
probable and are more unlikely to settle onto a configuration with well-separated kinks that 
can persist long oscillating transient. 
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FIG. 6. Likelihood of the observability of long transients starting from random initial 
conditions. The probability of observation of oscillating transients with at least 5, 10, and 50 
oscillations is obtained by numerical sampling of the space of initial conditions. The hypercube of 
initial conditions defined as ([0, m u ] x [0,p u ]) N has been sampled 10 4 times for each ring of length 
N. The likelihood of observing long oscillatory transients increases for longer rings. 

Figure [8] presents a summary of numerical simulations of a 16-gene ring that explore the 
characteristics of the two-kink traveling solution that make it observable. Firstly, random 
initial conditions settle with high likelihood onto configurations with two kinks separated 
by a distance k^ that is large enough to lead to a long-lived transient before annihilation. 
Figure |8^i shows the numerically observed proportion of initial conditions that settle onto 
configurations with different inter-kink distances. Note that 40% of the initial conditions 
lead to configurations with well-separated kinks (kd = 4, 5, 6) and the probability of reaching 
each of these distances is almost equal. 

As explained above, there is a relationship between the duration of the transient and 
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FIG. 7. Length of oscillating transients in rings of different sizes starting from compa- 
rable initial conditions. The same set of constants has been used for each ring. The insets show 
the initial conditions used corresponding to idealized 'up-down' solutions with 2 kinks separating 
domains Dl and D2. 



the inter-kink distance in the two-kink configuration. The further separated the kinks are, 
the longer the transient lasts, as shown in Figure [8^l, where we report the mean <T tr >oi 
the distribution of oscillating transient lengths T tr . Note the extremely large average of the 
transients ensuing from the symmetric (or near symmetric) two-kink arrangements with = 
5, 6. These long transients are responsible for the heavy tail of the probability distribution 
of the period lengths presented in Figure [HJo, which shows that the probability of observing 
a very long transient that oscillates for more than 550 periods is of comparable order to 
observing the stable 'up-down' fixed point directly. 

VI. COMPARISON TO OTHER UNIDIRECTIONALLY COUPLED 
RINGS: FLUX-GATE MAGNETOMETERS 

There are several systems with a directed ring topology where UPOs play a significant 
role in the dynamics. A recent example is the study of UPOs in relation to the emergence of 
chaotic behavior in rings of Duffing oscillators^. Closest to the dynamics of the gene system 
studied here are the rings of unidirectionally coupled magnetic cores 9 that have been not 
only studied numerically but also implemented experimentally by Bulsara et al. . In the 
antiferromagnetic case, the system of coupled magnetic cores is described by the following 
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(13%) - <T tt > = 78 periods 
(13%) - <TV> = 608 periods 
(14%) - <T tr > = 1587 periods 





H M I 


♦ 


k d = 




(MM 


♦ 


k d =\ 




Mill 


♦ 


k d = 2 


1 ! ! 1 ! ! 1 1 ! ! 


fill! 


♦ 


k d = 3 




H M 1 


♦ 


k d = 4 


t t i t i t 


H M 1 




k d = 5 


t ♦ ♦ t ♦ t ♦ t ♦ t 


'♦tit 


♦ 


k d = 6 



» 550 




50 100 150 200 250 300 350 400 450 500 550 
r tr in periods 



FIG. 8. Numerical analysis of oscillating transients starting from random initial condi- 
tions for a 16-gene ring, (a) The pie-chart represents the proportions of observed conformations 
collected from 10 4 initial conditions covering the hyper cube ([0,m u ] x [0,^])^. For each random 
initial condition, we measure the distance k d between the kinks after a burn-in time of five peri- 
ods. Combinatorially we expected that starting from random initial conditions each configuration 
should be reachable with comparable probabilities and this is confirmed numerically, as seen in the 
similar proportions of k d = 4, 5, 6. The proportions for the distances k d < 4 are pooled together 
due to the burn-in time of the simulation. The table also includes the mean of the distribution of 
numerically observed lengths of the oscillating transients for all the two-kink conformations for the 
16-gene ring, as represented below, (b) Distribution of the length of the transients showing a heavy 
tail due to the likelihood of observing the very long transients that ensue from the configurations 
with inter-kink distances k d = 4, 5, 6. 
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ODEtP: 

ij = —Xj + tanh (c (xj — 2^+1)) , (10) 

which falls under the framework of Eq. ^ . 

The bifurcations branching from the unstable symmetric solution follow from a Jacobian 
of the same form as that in Eq. ([7]). Hence the UPO characteristics follow the same pattern 
as those of the gene system. For even rings, the system shows symmetry breaking when 
c is increased giving rise first to the two degenerate dimerized solutions followed by the 
emergence of UPOs from the now unstable symmetric branch. For odd rings, there is also a 
similar cascade of Hopf bifurcations. The maximal number of Hopf bifurcations due to ring 
symmetry is (N/2 — 1) for even-numbered rings and (N — l)/2 for odd-numbered rings, as 
is the case for the generalized repressilator model. 

Such rings of magnetic cores have been implemented as real-world devices and the odd- 
numbered versions are used as sensitive magnetic sensors. Long even-numbered rings have 
also been built and their dynamics has been explored under noise perturbations and under 
random initial conditions^. It turns out that long oscillatory transients feature heavily in 
long even-numbered rings, even if the steady state is quiescentP^. 

VII. DISCUSSION 

The relevant dynamics of noisy systems with a limited horizon of observability is not 
always well defined by their long-term attractors as given by stability analysis. In such 
systems, long-lived transients can de facto dominate the observable behavior. Examples of 
such behavior are prevalent in a variety of engineering, physical and biological systems ^ 1 * 27 * 28 !. 
In this paper, we have analyzed the dynamics of the generalized repressilator model and have 
characterized the existence of long-lived oscillatory transients in even rings — systems that 
would otherwise be expected to behave as bistable switches based on strict stability analysis. 
These transients correspond to slowly decaying dynamics around particular unstable periodic 
orbits (UPOs) that emerge as part of a cascade of Hopf bifurcations. The UPOs have spatio- 
temporal symmetries that can be understood in terms of tilings and propagating kinks, 
corresponding to the moving walls separating different domains of the stable dimerized 
solution. In particular, the most probable long-lived observable transients correspond to 
propagation of kink pairs along the ring. We showed numerically that in the presence of such 
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UPOs the dynamics of the system is tunneled to its vicinity onto a constrained manifold from 
which it is difficult to escape, thus leading to observable very long-lived transient oscillations 
of hundreds of periods when starting from random initial conditions. The likelihood of such 
behavior can be understood from the fact that when a pair of kinks is well separated, they 
will propagate along the ring for a long time until they conflate to restore the fixed point 
solution. This process results in a heavy-tailed distribution of the duration of oscillatory 
transients. 

The gene rings analyzed here present analogies with other systems from classical physics, 
such as domain wall propagation in antiferromagnetic media or kink propagation in dimer- 
ized polymer chains, and is mathematically connected with classical models of nonlinear 
wave propagation in one-dimensional discrete lattices. However, most of previous models 
are based on undirected (bidirectional) coupling. Directed ring topologies with a coupling 
function that is bounded, monotonic and with a strong enough nonlinearity can lead to the 
emergence of UPOs with the right values of the coupling strength 4 19 . These UPOs can then 
potentially lead to observable long-lived oscillating transients. Recent studies have shown 
the importance of UPOs in similar directed ring systems, both as originators of transients 
in the route to chaos in rings of Duffing oscillators 7 and in the emergence of long-lived os- 
cillatory transients in experimental even-numbered rings of magnetic cores^. The linear 
stability analysis and the cascade of UPOs is similar in all those systems and follows from 
fundamental symmetry arguments^. However, the reachability and persistence of the tran- 
sient dynamics depends on the particular form of the coupling function and reflects global 
features of the flow in phase space that funnel the dynamics around the oscillatory behav- 
ior. Based on our interpretation in terms of kink propagation, it is likely that increasing the 
nonlinearity of the coupling function (while keeping it bounded and monotonic) will increase 
the likelihood of localized excitations. Hence, the existence of long-lived transients based 
on the propagation of kink pairs will become more likely. In the case of genetic circuits, 
this can be achieved by increasing the Hill coefficient h of the repressor. This can be done 
by the selection of particular gene promoters based on biochemical mechanisms with higher 
cooperativity. 

The design of topologically constrained systems is widely used in classical engineering and 
applied physics as an aid to controllable dynamics. The generic properties of such systems, 
which follow from fundamental results from group theory and dynamical systems^, provides 
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a good chance that the key underlying properties could survive under hardly predictable 
and classifiable biological environments. In this case, the symmetry of the ring arrangement 
leads to a strongly symmetric bifurcation analysis and also informs our understanding of the 
spatio-temporal solutions in terms of tilings and domain walls. However these topological 
considerations require that the elements in the ring be strictly identical. In real systems this 
condition is not fulfilled and the question remains if UPOs are dynamically observable when 
the individual elements are dissimilar, both in numerical simulations and in real-world built 
systems. The implementation of the ring of magnetic cores by Lindner and Bulsara 13 sug- 
gests that UPOs do feature in real-world applications. Our numerical analysis in Strelkowa 
and Barahona^ showed that UPOs are also observable in gene repressilator rings when 
the genes are not identical. They also feature heavily in the dynamics when the system is 
fundamentally stochastic, i.e., in the equivalent Markov process (Master equation). Such 
robustness allowed us to introduce the concept of a switchable synthetic genetic oscillator, 
operating under uncertainty and noise, based on the feedback operation around the UPOs 
associated with the long-lived transients. These examples and others indicate that the con- 
ceptual design based on simplified designs can still lead to a fundamental understanding, 
although system-specific characterizations in the presence of noise and uncertainty are essen- 
tial when applications to be considered. Such a perspective is of a potential value to provide 
ideas for transferable design unifying several disciplines, ranging from electrical engineering 
to synthetic biology. 
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Appendix A: Analytical conditions for the number of Hopf bifurcations 

Properties of the nonlinear coupling functions and the length of the ring N all deter- 
mine how many Hopf bifurcations and which UPOs are expected for a given parameter set 
{ci, C2, C3, C4} of the generalized repressilator model. 

Consider the case with the Hill repressor function (Eq. [2]). Then the characteristic poly- 
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nomial (J8I) becomes 



X 2 + (c 2 + c 4 )\ e + c 2 Ci (1 +aexp(ipi(N))) = (Al) 



where a = and £ = ^ with I = {0, .., N — 1}. 

The necessary condition for the Hopf bifurcation to occur is purely imaginary eigenvalues 
of the Jacobian. For a complex quadratic polynomial p(z) 

(a + iA)z 2 + 2(b + iB)z+ (c + iC) =0 (A2) 

the necessary condition for purely imaginary roots is^ 

C 2 - AbBC - Ab 2 c = 0. (A3) 

In our case, A = B = 0, b = (c2 + c&)/2, c = C2C4 (1 + a cos/3^) and C = c^c^a sin/3^, which 
leads to the condition: 

a 2 c 2 c 4: ( y l — cos 2 fa) — (c 2 + c 4 ) 2 (l + a cos fy) = 0. 

(A4) 

This equation determines which values of f3e(N) can satisfy the condition for Hopf bifurcation 
for a particular parameter set. 

Taking into account the symmetry of the cosine function and the fact that a>0, we 
recover the classical result!^ that the maximum number of solutions permitted by the ring 
symmetry is 

(N — l)/2 for odd rings 

(A5) 

(N/2 — 1) for even rings. 

However, depending on the parameter set some of those solutions might not be observed 
in numerical simulations. In our particular case, with Hill coefficient h = 2 and for the 
biologically relevant parameters^ c 2 = 0.12min~ 1 , C3 = 0.16 mm -1 , C4 = 0.06 mm -1 and 
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bifurcation parameter c\ G (0,3], we get that 



1 Cjv+i) 

2 2 



for odd rings 



(A6) 



[| (y — l)J for even rings. 



This set of solutions satisfies the bifurcation conditions for f3i(N) G (|7r, 2ir) and leads to 
the existence of UPOs in even rings with N > 6 as can be seen by using the bound 



< a 



hp. 



h+l 
m 



I 1 ) <h. 



Interestingly, for the same set of parameters {c;} but with a higher Hill coefficient h — 4, 



UPOs exist in smaller even rings and the criterion (A4) is fulfilled for N > 4. 
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